library(tidyverse)
library(WRS2)
library(cowplot)

#import dataset called "aid0018.xlsx"


#subset dataset appropriately
aid_eu15 <- subset(aid_0018, region == "EU15") #only western member states
aid_eu15_y14 <- subset(aid_eu15, year < 2015) #only years until 2014 for EU15
aid_ceec <- subset(aid_0018, region == "CEEC") #only CEECs
aid_ceec_member <- subset(aid_ceec, member == 1) #only years in which country is MS
aid_ceec_member_no_hr <- subset(aid_ceec, countrycode != "HR") #no Croatia
aid_ceec_member_y14 <- subset(aid_ceec_member, year < 2015) #only years until 2014 for CEECs
aid_eu26 <- subset(aid_0018, region != "Other") #take out those countries that are neither EU15 nor CEECs (Cyprus and Malta)
aid_eu26_no_hr <- subset(aid_eu26, countrycode != "HR") #no Croatia for regional subset
aid_eu26_y14 <- subset(aid_eu26, year < 2015) #only years until 2014 for regional subset
aid_eu26_member <- subset(aid_eu26, member == 1) #only years for which countries are MS for regional subset
aid_eu26_member_y14 <- subset(aid_eu26_member, year < 2015) #only years until 2014 for EU26


#create average aid over time by member state
aid_eu26_member <- aid_eu26_member %>%
  group_by(region, year) %>%
  mutate(mean_reg_aid = mean(totaid))

#robust t.test of differences in overall aid [footnote 34]
yuen(totaid ~ region, data = aid_eu26_member)

#robust t.test of differences in overall aid before 2014 [footnote 35]
yuen(totaid ~ region, data = aid_eu26_member_y14)

#Figure 1

#set boundaries for greyed-out part
xstart <- 2000
xend <- 2004

ggplot(aid_eu26, aes(x=year, y=mean_reg_aid, colour = region)) + 
  geom_line(size = 1.1)  + 
  geom_vline(aes(xintercept = 2004)) + 
  geom_rect(aes(xmin = -Inf, xmax = xend, ymin = -Inf, ymax = Inf),
            fill = "grey70", alpha = 0.009, colour = NA)  +
  scale_color_manual(values=c('black', 'grey20')) +  theme_light() +
  scale_y_continuous(expand=c(0.5,0)) +
  scale_x_continuous(expand=c(0,0)) + facet_wrap("region") +
  ylab("Average yearly aid as %GDP") + xlab("Year") +
  theme(legend.position = "none",
        plot.caption = element_text(size = 14),
        axis.title.x = element_text(size = 14),
        axis.title.y = element_text(size = 14),
        axis.text.x = element_text(size = 14),
        axis.text.y = element_text(size = 14),
        strip.text.x = element_text(size = 14)) +
  labs(caption = "Source: Personal elaboration on the EU Scoreboard data") 


#Figure 2

#make new df with only necessary elements
aid_eu26b <- data.frame(ID = aid_eu26_member$countrycode,
                        year = aid_eu26_member$year,
                        region = aid_eu26_member$region,
                        horaidperc = aid_eu26_member$horaidpc,
                        secaidperc = aid_eu26_member$secaidpc,
                        regaidperc = aid_eu26_member$regaidpc,
                        stringsAsFactors = T)

#generate mean aid by type over time by member state
mean_aid_eu26b <- aid_eu26b %>%
  group_by(ID) %>%
  mutate(mean_hor_aid = mean(horaidperc, na.rm = T)) %>%
  mutate(mean_sec_aid = mean(secaidperc, na.rm = T)) %>%
  mutate(mean_reg_aid = mean(regaidperc, na.rm = T))

#eliminate duplicated IDs
mean_aid_eu26b <- mean_aid_eu26b[!duplicated(mean_aid_eu26b[1]),]

#melt df so that only mean remains
aid_eu26b_melt <- data.table::melt(mean_aid_eu26b, id.vars = "ID")

#generate new df in which the means become percentages
aid_eu26b_melt2 <- aid_eu26b_melt %>%
  group_by(ID) %>%
  count(value, variable) %>%
  mutate(pct = value / sum(value) * 100,
         pct_label = scales::percent(pct)) 

# add ggplot
ggplot(data = aid_eu26b_melt2) + 
  geom_bar(aes(x = ID, y = pct, fill = variable),
           stat = "identity") +
  scale_fill_grey(labels = c("Horizontal", "Sectoral", "Regional"),
                  start = 0.7, end = .2) + theme_light() + 
  xlab("Country") + ylab ("Mean percentage of type of aid (2000-2018)") +
  guides(fill=guide_legend(title="Type of aid")) + 
  theme(plot.caption = element_text(size = 14),
        axis.title.x = element_text(size = 14),
        axis.title.y = element_text(size = 14),
        axis.text.x = element_text(size = 14),
        axis.text.y = element_text(size = 14),
        strip.text.x = element_text(size = 14),
        legend.text = element_text(size = 14),
        legend.title = element_text(size = 14)) +
  labs(caption = "Source: Personal elaboration on the EU Scoreboard data
       \nNote: Only years in which the states are EU members are considered")

#robust t.test in differences of type of aid (horizontal, sectoral, regional) [footnote 51]
yuen(horaidpc ~ region, aid_eu26_member)
yuen(secaidpc ~ region, aid_eu26_member)
yuen(regaidpc ~ region, aid_eu26_member)

#robust t.test in differences of type of aid until 2014 [footnote 51, but not shown in article]
yuen(horaidpc ~ region, aid_eu26_member_y14)
yuen(secaidpc ~ region, aid_eu26_member_y14)
yuen(regaidpc ~ region, aid_eu26_member_y14)


#Figure 3

#generate new df with mean over time by region
mean_aid_eu26b <- aid_eu26b %>%
  group_by(region, year) %>%
  mutate(mean_hor_aid_region = mean(horaidperc, na.rm = T)) %>%
  mutate(mean_sec_aid_region = mean(secaidperc, na.rm = T)) %>%
  mutate(mean_reg_aid_region = mean(regaidperc, na.rm = T))

#add ggplot
ggplot(mean_aid_eu26b) + 
  geom_line(aes(x=year, y=mean_hor_aid_region, col = "black"), size = 1.1) + 
  geom_line(aes(x=year, y=mean_sec_aid_region, col = "grey20"), size = 1.1) +
  geom_line(aes(x=year, y=mean_reg_aid_region, col = "grey70"), size = 1.1) +
  theme_light() +
  scale_x_continuous(expand=c(0,0)) + facet_wrap("region") + 
  scale_color_grey(name = "Type of aid", labels = c("Horizontal", "Sectoral", "Regional")) +
  ylab("Average yearly aid as % of total aid") + xlab("Year") +
  guides(fill=guide_legend(title="Type of aid")) +
  theme(plot.caption = element_text(size = 14),
        axis.title.x = element_text(size = 14),
        axis.title.y = element_text(size = 14),
        axis.text.x = element_text(size = 14),
        axis.text.y = element_text(size = 14),
        strip.text.x = element_text(size = 14),
        legend.text = element_text(size = 14),
        legend.title = element_text(size = 14)) +
  labs(caption = "Source: Personal elaboration on the EU Scoreboard data
       \nNote: Only years in which the states are EU members are considered")



#Figure 4

#generate new dataset for aid by objective 
aid_eu26c <- data.frame(country = aid_eu26_member$countrycode,
                        year = aid_eu26_member$year,
                        region = aid_eu26_member$region,
                        totaid = aid_eu26_member$totaid,
                        regdev = aid_eu26_member$regaid,
                        envaid = aid_eu26_member$env,
                        rdiaid = aid_eu26_member$rdi,
                        trainaid = aid_eu26_member$train,
                        sectaid = aid_eu26_member$sect,
                        smeaid = aid_eu26_member$sme,
                        cultaid = aid_eu26_member$cult,
                        empaid = aid_eu26_member$emp,
                        socaid = aid_eu26_member$soc,
                        heraid = aid_eu26_member$her,
                        rraid = aid_eu26_member$rr,
                        natdisaid = aid_eu26_member$natdis,
                        expaid = aid_eu26_member$exp,
                        closaid = aid_eu26_member$clos,
                        othaid = aid_eu26_member$oth,
                        stringsAsFactors = T)


#generate new df for means
mean_aid_eu26c <- aid_eu26c %>%
  group_by(region) %>%
  mutate(mean_regdev = mean(regdev, na.rm = T)) %>%
  mutate(mean_envaid = mean(envaid, na.rm = T)) %>%
  mutate(mean_rdiaid = mean(rdiaid, na.rm = T)) %>%
  mutate(mean_trainaid = mean(trainaid, na.rm = T)) %>%
  mutate(mean_sectaid = mean(sectaid, na.rm = T)) %>%
  mutate(mean_smeaid = mean(smeaid, na.rm = T)) %>%
  mutate(mean_cultaid = mean(cultaid, na.rm = T)) %>%
  mutate(mean_empaid = mean(empaid, na.rm = T)) %>%
  mutate(mean_socaid = mean(socaid, na.rm = T)) %>%
  mutate(mean_heraid = mean(heraid, na.rm = T)) %>%
  mutate(mean_rraid = mean(rraid, na.rm = T)) %>%
  mutate(mean_natdisaid = mean(natdisaid, na.rm = T)) %>%
  mutate(mean_expaid = mean(expaid, na.rm = T)) %>%
  mutate(mean_closaid = mean(closaid, na.rm = T)) %>%
  mutate(mean_othaid = mean(othaid, na.rm = T))


# melt data frame to long format
aid_eu26c_melt <- data.table::melt(mean_aid_eu26c[,c(3,19:33)], id.vars = "region")
aid_eu26c_melt1 <- aid_eu26c_melt[!duplicated(aid_eu26c_melt[2]),]

aid_eu26c_melt_eu15 <- subset(aid_eu26c_melt, region == "EU15")
aid_eu26c_melt_eu15 <- aid_eu26c_melt_eu15[!duplicated(aid_eu26c_melt_eu15[2]),]
aid_eu26c_melt_ceec <- subset(aid_eu26c_melt, region == "CEEC")
aid_eu26c_melt_ceec <- aid_eu26c_melt_ceec[!duplicated(aid_eu26c_melt_ceec[2]),]

#ggplot for EU15
type_aid_eu15  <- ggplot(data = aid_eu26c_melt_eu15) + 
  geom_bar(aes(x = variable, y = value),  fill =  "grey80", colour = "black",
           stat = "identity")  + theme_light() + coord_polar(clip = "off") +
  scale_x_discrete(labels = c("Regional \nDevelopment", "Environmental \nProtection", "R&D&I",
                              "Training", "Sectoral \nDevelopment", "SME", "Culture",
                              "Employment", "Social Support", "Heritage Conservation",
                              "R&R", "Compensation \nfor Natural Disasters",
                              "Export and \nInternationalisation", "Closure", "Other")) +
  theme(panel.border = element_blank(),
        axis.text.y=element_blank(),
        axis.ticks=element_blank(),
        legend.position="none",
        plot.title = element_text(hjust = 0.5),
        axis.text.x = element_text(angle=0, vjust = 0.5, hjust=1, size = 14),
        strip.text = element_text(size = 14)) +
  xlab("") + ylab ("") + facet_wrap("region")

#ggplot for CEECs
type_aid_ceec <- ggplot(data = aid_eu26c_melt_ceec) + 
  geom_bar(aes(x = variable, y = value, fill = variable), fill =  "grey60", colour = "black",
           stat = "identity")  + theme_light() + coord_polar(clip = "off") +
  scale_x_discrete(labels = c("Regional \nDevelopment", "Environmental \nProtection", "R&D&I",
                              "Training", "Sectoral \nDevelopment", "SME", "Culture",
                              "Employment", "Social Support", "Heritage Conservation",
                              "R&R", "Compensation \nfor Natural Disasters",
                              "Export and \nInternationalisation", "Closure", "Other")) +
  theme(panel.border = element_blank(),
        axis.text.y=element_blank(),
        axis.ticks=element_blank(),
        legend.position="none",
        plot.title = element_text(hjust = 0.5),
        axis.text.x = element_text(angle=0, vjust = 0.5, hjust=1, size = 14),
        strip.text = element_text(size = 14)) +
  xlab("") + ylab ("") + facet_wrap("region")

#create plot grid
type_of_aid <- plot_grid(type_aid_ceec, type_aid_eu15) 

#add captions
plot_type_of_aid <- ggdraw(add_sub(type_of_aid, 
                                   "Source: Personal elaboration on the EU Scoreboard data
                                 \nNote: Only years in which the states are EU members are considered",
                                   x = 1, hjust = 1.05, size = 14))


#Table 1 robust t.tests

#make aid by objective as perc of totaid
aid_eu26c2 <- aid_eu26c %>%
  group_by(region) %>%
  mutate(totaidpc = totaid/totaid*100) %>%
  mutate(regdevpc = regdev/totaid*100) %>%
  mutate(envaidpc = envaid/totaid*100) %>%
  mutate(rdiaidpc = rdiaid/totaid*100) %>%
  mutate(trainaidpc = trainaid/totaid*100) %>%
  mutate(sectaidpc = sectaid/totaid*100) %>%
  mutate(smeaidpc = smeaid/totaid*100) %>%
  mutate(cultaidpc = cultaid/totaid*100) %>%
  mutate(empaidpc = empaid/totaid*100) %>%
  mutate(socaidpc = socaid/totaid*100) %>%
  mutate(heraidpc = heraid/totaid*100) %>%
  mutate(rraidpc = rraid/totaid*100) %>%
  mutate(natdisaidpc = natdisaid/totaid*100) %>%
  mutate(expaidpc = expaid/totaid*100) %>%
  mutate(closaidpc = closaid/totaid*100) %>%
  mutate(othaidpc = othaid/totaid*100)

yuen(regdevpc ~ region, aid_eu26c2) #significant
yuen(envaidpc ~ region, aid_eu26c2) #ns
yuen(sectaidpc ~ region, aid_eu26c2) #ns
yuen(rdiaidpc ~ region, aid_eu26c2) #significant
yuen(empaidpc ~ region, aid_eu26c2) #ns
yuen(smeaidpc ~ region, aid_eu26c2) #significant


#Figure 5

#make df for mean by MS over time for CEECs only
aid_ceec <- aid_ceec %>%
  group_by(country, year) %>%
  mutate(mean_tot_aid = mean(totaid)) %>%
  mutate(mean_hor_aidpc = mean(horaidpc)) %>%
  mutate(mean_sec_aidpc = mean(secaidpc)) %>%
  mutate(mean_reg_aidpc = mean(regaidpc)) 

#add ggplot
ggplot(aid_ceec, aes(x=year, y=mean_tot_aid, group = countrycode)) + 
  geom_line(size = 1.1) + theme_light() +
  scale_y_continuous(expand=c(0.1,0)) +
  scale_x_continuous(expand=c(0,0)) + facet_wrap("country") +
  ylab("Average yearly aid as %GDP") + xlab("Year") +
  theme(legend.position = "none") +
  theme(plot.caption = element_text(size = 14),
        axis.title.x = element_text(size = 14),
        axis.title.y = element_text(size = 14),
        axis.text.x = element_text(size = 12),
        axis.text.y = element_text(size = 12),
        strip.text.x = element_text(size = 14),
        legend.text = element_text(size = 14),
        legend.title = element_text(size = 14)) +
  labs(caption = "Source: Personal elaboration on the EU Scoreboard data") 


#Figure 6

#add ggplot
ggplot(aid_ceec) + 
  geom_line(aes(x=year, y=mean_hor_aidpc, group = countrycode, col = "black"), size = 1.1) + 
  geom_line(aes(x=year, y=mean_sec_aidpc, group = countrycode, col = "grey20"), size = 1.1) +
  geom_line(aes(x=year, y=mean_reg_aidpc, group = countrycode, col = "grey70"), size = 1.1) +
  theme_light() +
  scale_y_continuous(expand=c(0.1,0)) +
  scale_x_continuous(expand=c(0,0)) + facet_wrap("country") +
  scale_color_grey(name = "Type of aid", labels = c("Horizontal", "Sectoral", "Regional")) +
  guides(fill=guide_legend(title="Type of aid")) +
  theme(plot.caption = element_text(size = 14),
        axis.title.x = element_text(size = 14),
        axis.title.y = element_text(size = 14),
        axis.text.x = element_text(size = 12),
        axis.text.y = element_text(size = 12),
        strip.text.x = element_text(size = 14),
        legend.text = element_text(size = 14),
        legend.title = element_text(size = 14)) +
  ylab("Average yearly type of aid over total aid") + xlab("Year") +
  labs(caption = "Source: Personal elaboration on the EU Scoreboard data")


#Robust ANOVA test of intra-bloc differences [footnote 65]
##totaid
t1way(totaid ~ country, aid_ceec_member)
lincon(totaid ~ country, aid_ceec_member) #for post-hoc comparisons
##horaid
t1way(horaidpc ~ country, aid_ceec_member)
lincon(horaidpc ~ country, aid_ceec_member) #for post-hoc comparisons
##secaid
t1way(secaidpc ~ country, aid_ceec_member)
lincon(secaidpc ~ country, aid_ceec_member) #for post-hoc comparisons
##regaid
t1way(regaidpc ~ country, aid_ceec_member)
lincon(regaidpc ~ country, aid_ceec_member) #for post-hoc comparisons


#Figure 7

#generate new df with median total and horizontal aid
comp_eu26_aid <- aid_eu26_member %>%
  group_by(country) %>%
  mutate(med_country_aid = median(totaid)) %>%
  mutate(med_country_horaid = median(horaidpc)) 


#add ggplot
ggplot(comp_eu26_aid) +
  geom_hline(yintercept = 65.57973, linetype = "dashed") + #median aid to be calculated separately
  geom_vline(xintercept = 0.475, linetype = "dashed") + #median horizontal share to be calculated separately
  geom_text(aes(x = med_country_aid, y = med_country_horaid, label = country),
            col = ifelse(comp_eu26_aid$region == "CEEC", "black", "grey70")) +
  theme_light() + 
  xlab("Median aid as %GDP (2000-2018)") + 
  theme(plot.caption = element_text(size = 14),
        axis.title.x = element_text(size = 14),
        axis.title.y = element_text(size = 14),
        panel.grid = element_blank()) +
  annotate(geom="label", x=1.4, y=90, 
           label="The dashed lines represent the EU-26 median \non the two dimensions",
           color="grey20", hjust = 1, size = 4.5) +
  ylab ("Median share of horizontal aid as % total aid (2000-2018)") +
  labs(caption = "Source: Personal elaboration on the EU Scoreboard data
       \nNote: Only years in which the states are EU members are considered") 


#Table 2


#pre- and post-accession FDI
fdi_ceec_04 <- subset(aid_ceec, year < 2004) #for all CEECs except BG and RO
fdi_ceec_07 <- subset(aid_ceec, year < 2007) # for BG and RO
fdi_ceec_18a <- subset(aid_ceec, year > 2003) #for all CEECs except BG and RO
fdi_ceec_18b <- subset(aid_ceec, year > 2006) # for BG and RO

#calculate mean FDI stock pre-accession for all CEECs except BG and RO for the 2000-2003 period
fdi_ceec_04 <- fdi_ceec_04 %>%
  group_by(country) %>%
  mutate(mean_fdi_pre = mean(fdi_inward_stock_tot)) 

#eliminate duplicates
fdi_ceec_04 <- fdi_ceec_04[!duplicated(fdi_ceec_04[1]),]
#only keep relevant columns
fdi_ceec_04 <- fdi_ceec_04[c(1,46)]

#calculate mean FDI stock pre-accession for BG and RO for the 2002-2006 period
fdi_ceec_07 <- fdi_ceec_07 %>%
  group_by(country) %>%
  mutate(mean_fdi_pre = mean(fdi_inward_stock_tot)) 

#eliminate duplicates
fdi_ceec_07 <- fdi_ceec_07[!duplicated(fdi_ceec_07[1]),]
#only keep relevant columns
fdi_ceec_07 <- fdi_ceec_07[c(1,46)]


#calculate mean FDI stock post-accession for all CEECs except BG and RO for the 2004-2018 period
fdi_ceec_18a <- fdi_ceec_18a %>%
  group_by(country) %>%
  mutate(mean_fdi_post = mean(fdi_inward_stock_tot)) 

#eliminate duplicates
fdi_ceec_18a <- fdi_ceec_18a[!duplicated(fdi_ceec_18a[1]),]
#only keep relevant columns
fdi_ceec_18a <- fdi_ceec_18a[c(1,46)]

#calculate mean FDI stock post-accession for BG and RO for the 2007-2018 period
fdi_ceec_18b <- fdi_ceec_18b %>%
  group_by(country) %>%
  mutate(mean_fdi_post = mean(fdi_inward_stock_tot)) 

#eliminate duplicates
fdi_ceec_18b <- fdi_ceec_18b[!duplicated(fdi_ceec_18b[1]),]
#only keep relevant columns
fdi_ceec_18b <- fdi_ceec_18b[c(1,46)]

